Terahertz conductivity of twisted bilayer graphene 
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Using terahertz time-domain spectroscopy, the real part of optical conductivity [(Ti(a;)] of twisted 
bilayer graphene was obtained at different temperatures (10 - 300 K) in the frequency range 0.3 - 
3 THz. On top of a Drude-like response, we see a strong peak in aiiui) at ~2.7 THz. We analyze the 
overall Drude-like response using a disorder-dependent (unitary scattering) model, then attribute 
the peak at 2.7 THz to an enhanced density of states at that energy, that is caused by the presence 
of a van Hove singularity arising from a commensurate twisting of the two graphene layers. 



Compared to single-layer graphene (SLG), where there 
are two non-equivalent lattice sites {A and B), bilayer 
graphene (BLG) has two SLGs stacked in the third di- 
rection. In the most common Bernal [AB] stacking of 
BLG, adjacent layers are rotated by 60°, so that the B 
atoms of layer 2 {B') sits directly on top of A atoms in 
layer 1 [A), and B and A' atoms are in the center of 
the hexagons of the opposing layer. Electrons can then 
hop between these two A sites with a hopping energy 
t^_. In the undoped case, though both SLG and BLG are 
gapless semi-metals, carriers in SLG exhibit linear dis- 
persion, while those in BLG show quadratic dispersion. 
An energy gap in SLG opens up due to finite geometry 
effects, but its control has proven to be unreliable [1]. On 
the other hand, the electronic gap in BLG can be reli- 
ably opened and controlled by an applied electric field, 
shown theoretically and demonstrated experimentally 

and promises interesting applications. Both SLG and 
BLG however, are sensitive to disorder. Hence, to realize 
graphene-based optoelectronic devices, an understanding 
of the temperature and disorder effects in the transport 
and spectroscopic properties of BLG is needed. Temper- 
ature and disorder-dependent conductivity of BLG have 
been derived theoretically [1, 6]. Experimentally, spec- 
troscopies (from terahertz (THz) to visible) and ultra- 
fast dyiramics of various flavors of graphene have been 
reported, such as SLG, few and many-layer graphene, 
and graphite (7l-[ll|. For example, Fourier-transform in- 
frared spectroscopy (FTIR) oir large-area SLG grown by 
chemical vapor deposition (CVD) revealed a Drude-like 
frequency dependence of the spectral density from THz 
to mid-infrared at different carrier concentrations [l^ . In 
addition, graphene plasmons, which lie in the THz range, 
are strongly coupled to the interband electronic transi- 
tions and decay by exciting interband electron-hole pairs 



[13|. Hence knowledge of graphene's electromagnetic re- 
sponse, as a function of disorder, in the THz frequency 
range is critical for applications such as graphene-based 
THz oscillators [li]. 

BLG grown by CVD also has a great tendency to twist. 
A typical 10 mm x 10 mm piece of CVD-BLG has been 
shown to be a collection of crystallites of twisted BLG 
with a distribution of different twisting angles (l5l . [3| ■ 
Twisting occurs when there is rotation between the top 
and bottom layers of BLG (see Fig. [T|). When there is 
rotation through a (twisting) angle 9 about an A {B') 
site in BLG, only a discrete set of commensurate airgles 



is allowed [17|: 
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where i = 0, 1, 2, .... Such rotation between graphene lay- 
ers have been observed as a Moire pat tern on graphite 
surfaces and recently in BLG |l5l|. Such twisting 
causes van Hove singularities (VHS) to develop near the 
Fermi energy, with the VHS energy scale being a strong 
function of 9, resulting in an enhanced density of states at 
those energies [l7j. Such enhancement in the density of 
states should show up in the conductivity spectrum. For 
example, for large twisting angles of 7.5°, 13.7° and 54.6°, 
anomalies in the real conductivity (Ti(cli) were seen in the 
visible region using contrast spectroscopy [19]. However, 
they have not been demonstrated in the THz regime. An 
accurate characterization of electrical and optical con- 
ductivities at THz frequencies of BLG, as a function of 
temperature, disorder and twisting, is thus needed, but 
has not been reported. 

In this Letter, we present THz time-domain spec- 
troscopy (THz-TDS) studies of twisted BLG at different 
temperatures (10 K - 300 K), to study its frequency- 
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FIG. 1: (Color online) Atomic arrangement of atoms in BLG, 
for twisting angle 6=21.8° (corresponding to i = 1). A [A') 
and B [B') are the sublattices of the first (second) layer. 
The black hexagon depicts the unit cell of the twisted BLG. 
Compare with Bernal (AB)-type stacking, where 61=60° (i.e. 
i = 0). 



dependent far-infrared conductivity. On top of a Drude- 
like response, we see a peak in the real conductivity. 
The overall Drude shape was analyzed using a disorder- 
dependent model, while the conductivity peak at 2.7 THz 
was attributed to an enhanced density of states at that 
energy, that is caused by the presence of low-energy VHS 
arising from a commensurate twisting of the graphene 
layers relative to each other. 

The samples studied here are large-scale BLG grown 
by CVD and deposited on z-cut quartz. Both contrast 
and absorption spectroscopies confirmed the sample to be 
a BLG film (2^. Our experimental set-up performs an 
average over the entire of the sample. Nevertheless, 
our data were able to discern the feature arising from 
Qi=2s ~ 1.161°, on top of a broad background produced 
by the disorder in the sample, as will be discussed later 
in this Letter. 

The transmission THz spectra of the BLG were mea- 
sured using a conventional THz-TDS system (TeraView 
Spectra 3000) with a Janis ST-IOO-FTIR cryostat. THz- 
TDS has proven to be a very useful noncontact tech- 
nique to study material properties such as dielectric re- 
sponse, complex conductivity and refractive index in the 
far-infrared range without the need for Kramers-Kronig 
analysis [2ll 22 1. The THz wave was generated and de- 
tected by photoconductive antennas fabricated on low- 
temperature-grown GaAs films. The aperture diameter 
is 7 mm, enabling accurate measurements of the low- 
frequency spectral components of the THz wave. For 
each sample or reference run, 900 THz traces were taken 
over 180 seconds. The sample holder was moved from 
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FIG. 2: (Color online) (a) Amplitude spectra of Sample 1, 
obtained from the Fourier Transform of the main pulse in the 
inset. Inset; THz-TDS signal of Sample 1 and reference at 
10 K. (b) Amplitude of complex transmission coefficient (or 
transmittance), at 10 K and 300 K. The solid line indicates 
the transmission amplitude when both sample and reference 
positions are vacuum. 



the sample to the reference position (and vice versa) by 
means of a vertical motorized stage with 2.5 /im reso- 
lution. The time-domain electric fields of a THz pulse 
transmitted through Sample 1 (Es{t) — BLG deposited 
on 1-mm thick z-cut quartz substrate from CrysTec, Ger- 
many), as well as through the reference {Er{t), bare z- 
cut quartz substrate) are shown in the inset of Fig. EJJa). 
Before BLG deposition, the substrates for sample and 
reference were carefully characterized by THz-TDS — 
their phase difference yields the thickness difference AL 
between the two substrates, which must be taken into 
account in our subsequent analysis (31. After the main 
pulse, a weaker pulse (etalon) appears due to multiple 
reflections in the z-cut quartz substrate. Since the main 
pulse and etalon pulse are well separated in the time do- 
main, we truncate the time-domain data to remove the 
etalon pulse. Subsequent data analysis was performed 
only on the main pulse without loss of validity. Fast 
Fourier Transform (FFT) was then performed to obtain 
the amplitude and phase at different spectral compo- 
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nents of the THz pulse. The FFT amplitude spectrum 
of the main pulse is shown in Fig. [IJa). The absorption 
of the THz pulse by the BLG is obvious, even though 
the sample is of atomic-scale thickness. Figure [2{b) 
shows the amplitude of the experimental transmission 
coefficient (or transmittance) T(aj), defined as the ra- 
tio between complex spectral field of the sample Es{<^) 
and reference Er{uj), for the BLG sample at 10 K and 
300 K. For both temperatures is almost frequency- 

independent with the value ~95 %. In the same figure is 
|T(a;)| when both sample and reference are vacuum — 
in this case deviates only 0.5% away from unity 

in the frequency range 0.3 - 3.0 THz, which will be the 
frequency window of our analysis. 

Theoretically, for a sample grown on a substrate, T(uj) 
can be written as 24l 
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2n{hsub + 1) cxp[iLud{n — l)/c] exp[—iu!AL{hsub — l)/c] 
(1 + n){n + fisub) + {n- l)insub 



h) ex'p[2iujdh/( 
_ (2) 

where h and TT-sub 8^re the complex refractive indices of 
BLG and z-cut quartz substrate, respectively, d {= I nm) 
is the thickness of the BLG [25], AL = —14 //m is 
the thickness difl:erence between sample and reference 
substrates (measured with a precision micrometer, and 
confirmed by THz-TDS data of the two bare substrates 
before BLG deposition), and c is the speed of light in 
vacuum. This expression takes into account the mul- 
tiple internal reflections inside the BLG sample, but 
does not include multiple reflections in the substrate — 
we need not take substrate reflections into account be- 
cause we have truncated the etalon pulse in our analy- 
sis. The complex refractive index ns„h of z-cut quartz 
was first measured with vacuum as reference at differ- 
ent temperatures, obtained to be fisub ~ 2.11 + 0.002i. 
This agrees with Ref. [i^ showing that z-cut quartz is a 
very good THz transparent material with a temperature- 
independent, and almost frequency-independent, refrac- 
tive index, in our frequency and temperature range. The 
complex refractive index h(uj) = n(u)) + ik(uj) is then 
extracted from Eq. ([2|) by numerical iteration, which is 
then used to calculate the complex optical conductiv- 
ity <7{u)) ~ cri(w) +ia2{i-o)^ where CTi(a;) — 2nkujeo and 
(T2(w) — {Soo — + fc^)a;£o7 being the free space per- 
mittivity, and high frequency dielectric constant Eoo — 8 
for graphene ^1]. However, the values of cr2(w) are very 
sensitive to the value of AL, due to the very small thick- 
ness of BLG compared to AL (~ fim). Hence we only 
discuss aiitj) in our subsequent analysis. Note that, for 
a very thin metallic film on an insulating substrate, the 
following assumptions can be used: n 3> Jisufc > 1 and 
dnoj/c ^ 1, and Eg. Q becomes the commonly- used 
thin-film expression [27| 

T{uj) = -I , - + "^"b exp[iujAL{iisub - l)/c] (3) 

where Zq is the free space impedance. The values of cri 



a 

O 

><, 




-20KFit 
-100 K Fit 
-300 K Fit tperp 



0^ — 
0.0 0.5 



1.0 1.5 2.0 

Frequency (THz) 



2.5 3.0 



FIG. 3: (Color online) Real conductivity ai{ij) of BLG of 
(a) Sample 1 and (b) Sample 2. Circles = data. Solid lines 
= simultaneous fits of the unitary scattering model to the 
10 K (black), 100 K (red) and 300 K (blue) data. Vertical 
axes on the right expresses ai{uj) in units of the minimum 
conductivity 8e^/(7r/i) as specified in Nilsson's model [l|. 



obtained from Eq. ^ are identical to that from Eq. 

Figure [Ita) shows cri(w) at 10 K, 100 K and 300 K 
of Sample 1. Note the very small difference between 
cri(a;,10 K) and tTi(a;,100 K). Superposed on top of 
a Drude-like response, is a strong peak centered at 
^2.7 THz. We first analyze the Drude-like background 
using a theoretical model developed by Nilsson et al. Q 
for unitary scatterers for Bernal BLG. The applicability 
of this theoretical model comes from the fact that disor- 
der broadens the low-energy features which would oth- 
erwise differentiate between the perfect Bernal-stacked 
BLGs and twisted BLGs. Therefore, we anticipate the 
robust validity of this model for the analysis of the Drude- 
like background. The model starts by considering a 
Hamiltonian of the BLG under the tight binding model. 
Within the T-matrix approximation for unitary scatter- 
ing, one derives the electron self-energy of BLG, which 
gives the Green's function in the presence of disorder via 
the Dyson equation. The conductivity is then calculated 
from the convolution of the Green's function elements 
(encoded in the kernal S), as a function of chemical po- 
tential fi, impurity concentration Ui, interlayer coupling 
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(hopping) t±_ and temperature, to be 



(Ti(a;) cx 



nh 



de 



(4) 

where np is the Fermi distribution function. The pre- 
factor 8e'^ /(irh) is the approximate minimal conductivity 
per BLG, whose exact value will depend on the actual 
distribution of impurities among the inequivalent sites 
of the A and B sublattices [ij]. In Fig. ^a), the 10 K, 
100 K and 300 K data were simultaneously fitted with the 
model via Eq. (g]), shown by black (10 K), red (100 K) 
and blue (300 K) solid lines. The resulting fitting pa- 



rameters were /i-^'* 



-0.044 eV, 



Jit 



0.00091, and 



t-^" = 0.049 eV. Note the small difference between the 
10 K and 100 K fits — consistent with data. In fact, 
the theoretical cri(w, 100 K) is smaller in magnitude than 
CTi(aj,10 K), showing that a spectral weight redistribu- 
tion has taken place. We were initially surprised that the 
100 K conductivity should be so similar to the 10 K con- 
ductivity, with the 300 K conductivity lying above both 
of them. These features, however, can be captured by 
the impurity-scattering model, lending strong credence 
to the validity of the model in explaining the THz data 
of Bernal BLG. In the above fittings, the fitted values 
of fi and Hi are consistent with Raman data of the same 
sample on the same substrate [2^ , with charged impurity 
concentration < lO^'^ cm^^ (corresponding to < 0.0026 
per BLG) and ^ = (-0.042 ± 0.012) eV, determined by 
the Raman peak positions of the G and 2D bands, and 
the intensity ratio of the 2D and G peaks f2^-[3l|. The 
unitary scatterer concentration of the sample is typically 
< 10~^, which was calculated using the ratio of the D 
and G peak intensities [32]. 

Figure ^h) shows the frequency dependence of the 
real conductivity, f7i(a;), at 20 K, 100 K and 300 K, 
of Sample 2. The simultaneous fits (solid lines) of the 



data to Eq. ([4]) now yield fitting parameters /i 



fit 



-0.012 eV, * = 0.00071, and t^" = 0.076 eV. Once 
again the fitted is consistent with Raman data, where 
fi = (—0.012 ± 0.008) eV, and similar impurity concen 
tration as Sample 1 . Note that both fitted values of t^ 
are smaller than the value t±=0.27 eV for Bernal BLG. 
In fact, for a single monodomain of twisted BLG, the 
interlayer hopping is angle dependent, but for small an- 
gles it can be approximated as t^j_ OAt± ~ 0.1 eV (s^ . 
and <0.1 eV for larger 6''s (larger 9 implies a larger 
separation between the layers, hence smaller interlayer 
hopping). Hence the value of obtained (w50-70 meV) 
could be an average of interlayer couplings from all possi- 
ble twisting angles in the sample, re-expressed in the form 
of perfect Bernal stacking. The fit to a theory based on 
an "effective" Bernal BLG only works because disorder 
broadens all the features which would otherwise distin- 
guish Bernal from twisted BLG at low energies, namely, 
the presence of VHS in the density of states of twisted 
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FIG. 4: Schematic dispersion of a (a) Bernal BLG and (b) 
twisted BLG. The shaded region indicates the states that have 
been broadened by disorder in the sample, hence our data 
for twisted BLG could be fitted with a theory developed for 
Bernal BLG, but with a smaller value of t±. 



BLGs [17] (see Fig. S]). The difference in absolute val- 
ues of (Ti (w) between samples 1 and 2 is consistent with 
sample-to-sample variations observed in other works on 

BLG [ail. 

A close inspection of ai{uj) reveals the presence of a 
peak that appear on top of the background signal. In 
twisted BLG, VHS develop near the Fermi energy, which 
results in an enhanced density of states fl3|- The energy 
scale of such VHS depends sensitively on the twisting 
angle 9, given by 



Ev. 



hs 



3a 



sin(e'/2)| - 2t^ 



(5) 



where vp = 1.0 x 10^ m/s [33 is the Fermi velocity and 
a«2.46 A the lattice constant. We observed a strong peak 
at ~2.64 THz, whose presence is reproducible from sam- 
ple to sample. Its position is consistent with the second 
non-zero E^hs, computed from Eq. ([5]) to be 2.77 THz, 
and corresponds to 6*28 — 1.161° (from Eq. |T])). Note 
that the first non-zero Eyhs of 0.89 THz, arising from 
6*29 = 1.121°, is not visible from Fig. [31 Theoretical den- 
sity of states calculations [33] show that, for 027=1-20°, 
the van Hove peaks are still barely visible, whereas for 
6*30=1.08° the VHS have disappeared [s^. Also, since 
disorder builds up the density of states near the Dirac 
point, the VHS are broadened by being in the middle of 
a continuum of disordered states. These factors may ex- 
plain our inability to see any clear feature near 1 THz. 
Note that our 2.7-THz peak is robust against the type of 
windowing function we used before performing FFT. Be- 
sides the conventional windowing functions, we also con- 
structed an asymmetric windowing function that is tai- 
lored to the shape of our asymmetric time-domain wave- 
forms (36| — all yielded the 2.7-THz peak. This 2.7-THz 
(^11 meV) Eyfis is also consistent with scanning tunnel- 
ing microscopy and spectroscopy (STM/STS) works on 
twisted BLG [1^,113. 

The Raman data on the same samples also gave infor- 
mation about the twisting [2^. The position of the G 
peak shows the samples to be slightly p-doped [sij . The 
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blueshift of the 2D peak implies the existence of twist- 
ing, and the value of the 2D peak width implies a twist- 
ing angle < 5° The consistency of these results 
across different positions of the samples implies a well 
defined twisting angle in our samples. Hence our THz 
data, besides being consistent with Raman data in the 
same samples, points out the exact twisting angle, and 
shows the effect on twisting on the optical conductivity. 

In conclusion, we have studied the far-infrared dielec- 
tric response of bilayer graphene at different tempera- 
tures by THz-TDS. On top of a Drude-like response, we 
observed a peak in the real part of optical conductivity. 
The overall Drude shape was analyzed using a disorder- 
dependent model, while the conductivity peak at 2.7 THz 
was attributed to an enhanced density of states at that 
energy, that is caused by the presence of a low-energy van 
Hove singularity arising from a commensurate twisting of 
the top graphene layer relative to the bottom layer. A 
unified theory that considers the effect of both disorder 
and twisting on BLG conductivity is clearly desired. 
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